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We construct a simple two-phase equation of state intended to resemble that of compressed 
baryon-rich matter and then introduce a gradient term in the compressional energy density to 
take account of fintie-range effects in non-uniform configurations. With this model we study the 
interface between the two coexisting phases and obtain estimates for the associated interface tension. 
Subsequently, we incorporate the finite-range equation of state into ideal or viscous fluid dynamics 
and derive the collective dispersion relation for the mechanically unstable modes of bulk matter in 
the spinodal region of the thermodynamic phase diagram. Combining these results with time scales 
extracted from existing dynamical transport simulations, we discuss the prospects for spinodal phase 
separation to occur in nuclear collisions. We argue that these can be optimized by a careful tuning 
of the collision energy to maximize the time spent by the bulk of the system inside the mechanically 
unstable spinodal region of the phase diagram. Our specific numerical estimates suggest cautious 
optimism that this phenomenon may in fact occur, though a full dynamical simulation is needed for 
a detailed assessment. 

PACS numbers: 25.75.-q, 81.30.Dz, 64.75. Gh, 64.60.an 



I. INTRODUCTION 

The phase structure of strongly interacting matter 
presents a focal point for current theoretical and ex- 
perimental investigations. In particular, the Relativis- 
tic Heavy Ion Collider (RHIC) at BNL is preparing for 
a beam energy scan that aims to identify signals of the 
expected critical point and the CBM experiment at the 
future Facility for Antiproton and Ion Research (FAIR) 
at GSI will explore the properties of compressed bary- 
onic matter and search for the expected first-order phase 
transition. 

On the theoretical side, the situation is far from clear. 
Whereas lattice QCD calculations [J 0] find that the 
deconfinement phase transformation is of the crossover 
type at vanishing chemical potential, (i — 0, they have 
inherent difficulties treating finite /x values and any pre- 
dictions in the baryon-rich domain arc still very uncertain 
@; S Q • Even the very existence of a critical point has 
recently been called into doubt [fj. Experimental infor- 
mation would therefore be invaluable. 

However, it will be no easy task to extract the ther- 
modynamic phase structure from nuclear collision exper- 
iments. In addition to the inherent problems arising from 
the smallness of the collision system (which renders its 
spatial configuration far from uniform) and its rapid evo- 
lution (which prevents global equilibrium from being es- 
tablished) , the experimentalist is faced with the problem 
that there exists yet no suitable dynamical model with 
which to simulate the collisions for the purpose of anti- 
cipating the observable effects of the phase structure. 

This crucial point deserves elaboration: First of all, 
the basic theory, quantum chromodynamics, is currently 
tractable only in cither the perturbative limit of hard 
elementary processes or in the thermodynamic limit at 
vanishing (or small) net baryon density. Any dynamical 



transport treatment of nuclear collisions must therefore 
involve a considerable degree of modeling. 

Ideally, one would devise a transport model that ex- 
plicitly treats the dynamics of the microscopic degrees of 
freedom in the system, which change from being partonic 
in the deconfined sector to being hadronic in the confined 
sector. Unfortunately, it has yet not been possible to de- 
velop such a description, even for static scenarios. Nev- 
ertheless, a variety of microscopic transport models have 
achieved considerable success with regard to calculating 
(and reproducing) observables for high-energy collisions 
over a large range of energies. However, their thermo- 
dynamic properties are (yet) inadequate. For one thing, 
they usually lack detailed balance (as is often well justi- 
fied in the context of the dynamical processes for which 
they are intended) and therefore they are inadequate 
for thermal equilibrium. Furthermore, these models do 
not (as of yet) incorporate a first-order phase transition. 
Therefore, at this point, they appear to be unsuitable for 
simulations that aim to bring out the dynamical effects 
of a phase transition and elucidate their observability. 

Considerable success has been obtained as well with 
macroscopic models within the framework of fluid dy- 
namics. These models have the practical advantage that 
the underlying microscopic degrees of freedom do not en- 
ter explicitly, the state of the system being described 
merely through its local energy and charge densities (and 
the associated currents) with the interactions entering via 
the equation of state and the transport coefficients. [Of 
course, in order to make contact with experiment, such a 
treatment must ultimately convert the macroscopic infor- 
mation into hadrons by a suitable freeze-out prescription, 
but this occurs at densities well below the phase transi- 
tion region and is well developed.] Thus fluid dynamics, 
especially ideal fluid dynamics for which the transport 
coefficients vanish, posesses a very close relationship be- 
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tween the dynamics and the underlying object of study, 
the equation of state. 

However, a closer analysis reveals that standard fluid 
dynamics has certain inherent problems in the presence of 
a first-order phase transition. Of particular importance 
is the fact that standard fluid dynamics is strictly local 
which leads to both static and dynamic shortcomings, as 
we shall now discuss. 

With regard to the former, imagine that two thermody- 
namically coexisting bulk systems are brought into con- 
tact along a common interface. In a realistic description, 
a diffuse interface would develop, with the various densi- 
ties changing smoothly from one bulk value to the other, 
and there would be an interface tension. By contrast, 
when the equation of state is strictly local, the inter- 
face will be sharp (so the various densities will change 
abruptly from one bulk value to the other across the in- 
terface) and there will be no interface tension. Such a de- 
scription would not be adequate for finite systems such 
as blobs of matter produced in a collision, whose sizes 
are determined primarily by the surface tension and for 
which much of the matter is located in the diffuse surface 
region. While this generic shortcoming may be less seri- 
ous for high-energy collisions, where the matter is being 
torn apart due to the rapid longitudinal stretching, it is 
expected to play a significant role at the lower collision 
energies relevant for the exploration of the deconfinement 
phase transition. 

As for the dynamics, consider the evolution of nearly 
uniform matter that has been prepared in a state of ex- 
pansion at a density just above the phase coexistence re- 
gion. The system would then continue its expansion and 
the associated phase point would soon enter the phase 
coexistence region in which uniform matter is thermo- 
dynamically metastable. While this would pose no par- 
ticular problem as long as the deviations from unifor- 
mity remain small, the further expansion would drive the 
phase point into the region of spinodal instability, where 
uniform matter is both thermodynamically and mechan- 
ically unstable (the speed of sound is imaginary). As a 
result, density undulations would (and should) become 
amplified. 

This scenario is familiar from many areas of physics 
and it has been studied both theoretically and experi- 
mentally for a variety of substances [7|, [8( . Generally, the 
associated collective dispersion relation (which in this sit- 
uation gives the growth rate as a function of the wave 
number, 7^) exhibits a maximum, thus leading to prefer- 
ential amplification of certain modes and the appearance 
of a characteristic length scale in the ensuing phase sepa- 
ration. This remarkable phenomenon, known as spinodal 
phase decomposition, is an indicator of a first-order phase 
transition. It was found to present a powerful means 
for the experimental exploration of the nuclear liquid- 
gas phase transition d, Q , because the unstable dilute 
bulk matter tends to condense into fragments of simi- 
lar sizes, a highly non-statistical outcome that is easy to 
identify in the event analysis. This success has given rise 



to the hope that spinodal decomposition could be useful 
as well for probing the confinement phase transition and 
some explorations of possible experimental signals have 
already been made [jjl El El El • 

Naturally, since standard fluid dynamics is local, so 
is its collective dispersion relation, u>k = v$k. Conse- 
quently, inside the spindoal region of the phase diagram, 
the growth rate will increase monotonically with the wave 
number of the undulation. Thus jk will not display a 
maximum and the characteristic spinodal decomposition 
phenomenon would not be develop, as density irregulari- 
ties of ever smaller scale would be amplified at ever larger 
rates. In ideal fluid dynamics, this problem would be 
computationally intractable (and in fact mathematically 
meaningless) |14j . The inclusion of viscosity would mod- 
ulate the dispersion relation and cause the growth rate to 
approach a constant value for large k. While this would 
facilitate the numerics, the monotonic growth of jk would 
still preclude the occurrence of a spinodal decomposition. 

Because of the considerable potential for fluid dynam- 
ics as a tool for obtaining insight into the phase transition 
dynamics, we address here this generic shortcoming. Our 
main purpose is to illustrate the advantages of remedying 
this problem and, along the way, make rough estimates 
for various key quantities. Since the problem arises from 
the local nature of fluid dynamics, we introduce a finite 
range into the treatment by means of a gradient term in 
the compressional energy density. The specific model de- 
veloped here is intended to serve mainly as a framework 
for illustrating the effect of incorporating a finite range 
into the dynamcial desccription and the specific quantita- 
tive results should be regarded as correspondingly rough. 

We first (Sect. |TT| construct a somewhat schematic 
equation of state for uniform matter, trying to incor- 
porate the most essential features expected; it should 
be considered as merely a temporary substitute sub- 
ject to refinement. Subsequently (Sect. IIII[) the finite 
range is introduced by means of a simple gradient term. 
Then (Sect. [IV]) we consider the equilibrium interface be- 
tween two bulk systems, a property that could not be 
addressed with the standard treatment due to its strict 
locality, and obtain expressions for the associated inter- 
face tension. The collective modes in bulk matter are 
then treated (Sect. |V| and we consider particularly the 
spinodal growth rates which now display the characteris- 
tic features known from other substances. Finally (Sect. 
IVip , on this basis, we discuss the phase transition dynam- 
ics expected for the planned nuclear collision experiments 
and the prospects for spinodal decomposition to actually 
occur. 



II. BULK MATTER EQUATION OF STATE 

We wish to employ an equation of state that is suitable 
for numerical illustrations. For this purpose, we design 
a schematic model that is a generalization of a classi- 
cal gas in a density-dependent mean field. The resulting 
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equation of state has certain generic deficiencies and the 
results should therefore not be taken at face value. [For 
example, there are no bosonic degrees of freedom.] It 
would of course be of interest to repeat the present anal- 
ysis as more refined descriptions become available. 

The equation of state provides the thermodynamic 
properties of bulk matter, i.e. uniform matter of suffi- 
cient spatial extension to render finite-size effects (includ- 
ing those from any surfaces) insignificant. In the micro- 
canonical representation, the state of the system is given 
in terms of the basic mechanical densities, the (baryon) 
charge density p and the energy density e which we take 
as a thermal energy, n = \dpT, plus a compressional en- 
ergy, wq(p), where d is an adjustable parameter (equal 
to three for an ordinary gas). Further adjustable param- 
eters appear in the compressional energy density, wq (p) , 
which is specified in App. [SJ 

The key thermodynamic quantity is the entropy den- 
sity a(e, p), which we express in terms of the entropy den- 
sity <t(k, p) for a generalized ideal classical gas of density 
p and thermal density k, 
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a(s,p) = a(e -w (p),p) 
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where we have defiend the thermal density as 



Pr{e,p) 
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with T{e, p) — |[e — w (p)]/p (see below). The Lagrange 
coefficients (3 — l/T and a — —p/T are given by 



(3(e,p) = d £ a(e,p) = d K a(n - w (p), p) 
a(s, p) = d p a(e, p) = d p a(n - w , p) = c 



(3) 
(4) 



with a K = 5 k ct(k,p), a p = d p a(n,p), and w' = d p w (p). 
Thus, the temperature and the chemical potential are 
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p(e, p) 



2 e-wp{p) 
d p 



2 K 
d p ' 



(5) 



-aT 



UTln-?-+w' (p) . (6) 
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The pressure and the enthalpy density may be obtained 
subsequently, 



p(e, p) = aT - e + pp = \dpT - w + pw' , (7) 
h(s, p) = p + e = \dpT + pw' . (8) 



Two bulk systems with the densities (si, Pi) and 
{s2, pi) are in mutual thermodynamic equilibrium iff the 
total entropy is stationary under arbitrary exchanges of 
energy and charge, yielding the requirement that they 
have equal temperatures, chemical potentials, and pres- 
sures: (3i = 02 = Po, oil = a 2 = a , pi = p 2 = Po- Thus 
phase coexistence requires that the gradient of o{e,p), 
(<7 e , dp) = (d £ <7, d p a), be the same at the two phase points 
and, furthermore (since p = T[a— (3e — ap]), that the tan- 
gent to <r(e, p) at these two points be common. Further- 
more, local thermodynamic stability requires that the 
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FIG. 1: The equation of state pt{p). The pressure p as a 
function of the density p for a range of temperatures, T/T c = 
0, 4, |, |, 1, | , |, obtained with the adopted model. The 
phase coexistence (solid) and the spinodal (dashes) bound- 
aries are indicated; they coincide at the critical point (dot). 



second variation of the entropy be positive under such 
exchanges, yielding the requirement that the curvature 
matrix of tr(£, p) be positive definite. Consequently, the 
region of spinodal instability is delineated by the occur- 
rence of a vanishing curvature eigenvalue. 

In the canonical representation e is replaced by T, and 
the free energy density is then of special interest, 



f T {p) = s T (p) -Ta T (p) = pr(p)p- Pt(p) 

= pTln pT + w Q (p) , 

Pt 



(9) 



where the subscript T indicates that the quantity is ob- 
tained at the specified temperature. We also note that 
the slope of the free energy density is the chemical po- 
tential, d p f T {p) = Pt(p)- 

For wq(p) we use an interpolated form that produces 
a first-order phase transition (see App. [A}. To achieve 
a semi-quantative correspondence with expectations for 
actual baryon-rich matter, we adjust the parameters such 
that the coexistence densities at T = are p\ = 3po and 
P2 = 8po, where po « 0.153 fm~ 3 is the nuclear satu- 
ration density; the associated zero-temperature specific 
heat is then w (p 2 ) — w (pi) — 590 MeV/fm 3 . Further- 
more, the value d — 5.5 yields a critical temperature of 
T c = 170 MeV; the critical density is then p c = 4.70 pq. 
Other valus of particular interest are listed in Table [U 
(These values are of course somewhat arbitrary but will 
serve well for illustrative purposes.) The resulting equa- 
tion of state, pt(p), is shown in Fig. [fl while Fig. [5] dis- 
plays the associated phase diagram expressed in terms 
of the mechanical densities p and e. The more famil- 
iar (p, T) phase diagram, for which the energy density 
e has been replaced by the temperature T, is shown in 
Fig. [3] It is important to recognize that whereas the 
transformation from e to T is always unique, the reverse 
transformation is triple- valued in the presence of a phase 
transition: Any (p, T) phase point inside the phase coex- 
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Compression p/p s 

FIG. 2: The phase diagram in the p-e plane, as obtained 
from the equation of state (Fig.rj}, with the phase coexistence 
boundary (solid, red), the (isothermal) spinodal boundary 
(short dashes, green), and the critical point (dot) indicated. 
The hadronic freezeout line (lower left) is included for refer- 
ence (from Ref. [lH). Also shown are the two functions wh(p) 
(dashed blue curve) and wq(p) (dashed red curve) between 
which the compressional energy wo(p) (solid curve) is inter- 
polated, as well as three isentropic phase trajectories (dot- 
dashed curves), for which p8e — (e+p)Sp. 



istence region in Fig. [3] could arise any of three different 
(p, e) phase points in Fig. [2] 

At a given temperature T, bulk matter at the two 
different densities p\ and p2 are in mutual thermo- 
dynamic equilibrium if the corresponding tangents of 
/t(p) are common: the two chemical potentials are then 
equal since pr{p) = 5 p /t(p), and the relation pr{p) = 
Pt(p)p — Jt(p) guarantees that also the two pressures 
match. Thus phase coexistence at T = requires that 
the tangents of wo(p) at the two densities be common. 
[We have used this propertiy to guide our choice of mean 
field.] As the temperature is increased, the difference be- 
tween the two coexistence densities will steadily shrink 
until they coincide at the critical temperature T c . 

At supercritical temperatures (T > T c ) the pressure 
increases steadily with T, d p pr > 0, whereas its be- 
havior is undulatory at subcritical temperatures: when 
the density is increased from the lower to the higher 
coexistence density, the pressure exhibits first a maxi- 
mum and then a minimum. The associated densities pa 
and pb where pt(p) is stationary delineate the region 
of mechanical instability, within which d p px(p) is nega- 
tive. Since d p pT(p) = f T + PWq , the spinodal boundary 
densities at T = are determined by pwq(p) = and 
w'q(p) is negative in between. The region of mechanical 
instability shrinks steadily as T is increased and disap- 
pears at T c , which is thus determined by the condition 
|T C + Pm'WQ(pm) = 0, where p m is the density at which 
PWq(p) is most negative. The adopted compressional en- 
ergy density wo(p) is depicted in Fig. [1] together with the 
coexistence and spinodal boundaries. 

Dynamical transport calculations suggest that the ex- 




Compression p/p s 



FIG. 3: The phase diagram in the p-T plane indicating the 
phase coexistence boundary (solid), the isothermal spinodal 
boundary (long dashes), the isentropic spinodal boundary 
(short dashes), and the critical point (dot). Also shown are 
several isentropic phase trajectories, for which p5e= (e+p)Sp. 

pansion stage in a nucleus-nucleus collision proceeds in 
an approximately isentropic manner [16j| . i.e. the en- 
tropy per (net) baryon remains nearly constant. Since 
T p 2 5{a I ' p) = (pSe — ppSp) — (h5p — ppSp) — pSe — hSp, 
the isentropic trajectories in the (p, e) phase plane are 
characterized by p6e = hSp. Figs. and [3] display several 
such isentropic phase trajectories and they are seen to not 
be noticeably affected by the presence of the phase tran- 
sition. This feature brings out the fact that the locations 
of the boundaries for thermodynamic and mechanical in- 
stability, including the critical point itself, result from a 
rather subtle interplay between the underlying interac- 
tions. One may therefore expect that the overall phase 
evolution obtained in a dynamical transport calculation 
is not very sensitive to the specific phase structure. 



III. GRADIENT CORRECTIONS 

The above thermodynamics discussion applies to bulk 
matter, i.e. large and uniform systems. In heavy-ion 
physics, the systems encountered are neither and it is 
therefore practically important to extend the treatment 
to systems whose densities vary with the location, s(r) 
and p(r), where we use a tilde over a quantity as a re- 
minder that it pertains to a non-uniform system. 

As a simple way to take approximate account of finite 
range effects, we employ a gradient correction in the com- 
pressional energy. (A gradient term was also employed in 
recent hydrodynamical studies of the hadron-quark first- 
order phase transition [l7j|.) Accordingly, we write the 
local interaction-energy density on the form 

tB(r) = w (p(r)) + iC(V/5(r)) 2 (10) 

= Mp(r)) + WeA^-) 2 • 
\ Pg / 
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It is convenient to write the strength of the gradient term 
on the form C = a 2 e g /p 2 , where p g is a characteristic 
charge density and e g is a characteristic energy density. 
Since we are here particularly interested in the dynam- 
ics in the phase transition region, we choose the phase 
point (p g ,£g) to be in the middle of the phase coexis- 
tence region, p g = p c = 4.70 p and e g = e T = Tc/2 {Pc) = 
\{w {p c )+e c ) = 561MeV/fm 3 . The strength of the gra- 
dient term is then governed by the length a which we 
consider to be somewhat adjustable. Our present calcu- 
lations have been made with a = 0.2 fm. 

The introduction of the gradient term leads to gradient 
corrections in the expressons for the various thermody- 
namic quantities. In order to derive those, we start from 
the entropy density, which we assume to still have the 
form <r(r) = a(k(r), p{r)) = cr(e(r) — w(r), p(r)), where 
k(r) = e(r) — w(r) is the local thermal energy density. 
A variation of the total entropy S[e(r), p(r)] = Jdr a(r) 
then yields the local Lagrange coefficients (3 and a, 



p{r) 
a(r) 



5S 
5e{r) 

6S 
8p{r) 



= K(n(r),p(r)) = l/f(r) , (11) 
= °a p (K(r),p(r))-p(r)w' (p(r)) 



CV0(r) ■ Vp(r)) = -£(r)/T (r) . (12) 



Using that the entropy density gradient is then given by 
Vcr(r) = P(r)Ve(r)+a(r)Vp(r))-CV(P(r)(Vp(r)) 2 ), 

we see that the following expression for the local pressure, 
p(r) = &(r)f(r)-e(r) + p(r)p(r) + C{Vp(r)) 2 , (14) 
leads to the relation 



T(r) 



-e(r)V/3(r) - p(r)Va(r) , 



(15) 



which can be regarded as a generalization of the familiar 
thermodynamic relation 8{p/T) = —eS(3 — pSa. This 
relation ensures that p(r) will be constant whenever T(r) 
and p,(r) are. We also note that the gradient correction 
to the compressional energy migrates directly into the 
free energy density, 



f T (r) = K T (p(r))+w(r)-T?r(K T (p(r)),p(r)) , 



= Mp(r)) + iC(Vp(r)) 2 . 



IV. INTERFACE EQUILIBRIUM 



(16) 



Once the finite-range effects have been included in the 
thermodynamics, one may treat the interface between 
two coexisting phases. For this purpose, we consider a 
semi-infinite geometry with the two coexisting systems 
having a planar interface perpendicular to the x direc- 
tion. The coexistence values of temperature, chemical 
potential, and pressure are denoted by Tq, po, and pq. 
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FIG. 4: The specific interface tension -y Tg as a function of the 
coexistence temperature To for various values of the range a. 



We first note that global equilibrium, including equi- 
librium between two bulk systems with a common inter- 
face, requires that the total entropy S be constant under 
variations Se(x) and Sp(x) that conserve the total energy 
E = f dxe(x) and the total (net) charge B = J dx p(x), 



= SS - (3 SE - a a 5B 

= 5 dx [a(x) — Poi{x) — aop(x)] 



(17) 



= dx 



{[/3(x) - p ]Se{x) + [a(x) - a }6p{x)} 



thus implying spatial constancy of the temperature and 
the chemical potential, (3(x) = (3q and a(x) = ao, as one 
should expect. 

Therefore, assuming that the temperature is constant, 
(3(x) = (3q = 1/To, it is convenient to work in the 
canonical framework and our analysis is then similar 
to that carried out by Ravenhall et al. [18]. With the 
temperature given, the local density p(x) determines 
the local energy density, e(x) — |p(x)T + w(x), and 
the local entropy density is then in turn determined, 
a(x) = a(e(x) — w(x) : p(x)) = a(^p(x)T , p{x)). The 
local free energy density is then readily obtained, 

/>) = i(x)-T a(x) = fTMx)) + \C{d x p{x)) 2 , (18) 

where /t(p) is the free energy density in bulk matter at 
temperature T and density p (see Sect. [TTJ) . The cor- 
responding bulk chemical potential is pt{p) = dp/rip), 
while the bulk pressure is pr{p) — Pt(p)p — It{p)- 
The condition for equilibrium can now be expressed as 



= 5 dx 



f(x) - p p(x) 



(19) 



= J dx [p Tn (p(x)) - Cd 2 x p(x) - no] Sp(x) 



which then requires 

Cd 2 x p{x) = p T Mx))-M = d P Af(p(x)) . (20) 
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FIG. 5: Top: The surface profile p(x) for T = (solid) and 
T = |T c (dashed), using as a reference for x the position xq 
where the chemical potential equals the coexistence value of 
the bulk chemical potential, pr(p(x = Xq)) = po(T). The 
limiting (coexistence) densities are shown by the horizonthal 
lines, while the bottom curve is the interface location function 
g(x) for T = 0. Bottom: The mean location x of the interface 
and its width b as functions of temperature. Also shown is 
the profile skewness parameter 73 = ((x — x) 3 )/b 3 . 



Here A/(p) is the difference between the free energy den- 
sity of a uniform system of density p, /t (/°)j an d the 
corresponding "Maxwell" free energy density, defined as 
the free energy density along the common tangent, 



JtM = frM + nip-Pi) < h {p) 



(21) 



where pi refers to either one of the two coexistence den- 
sities. Thus Af(p) can be thought of as the free energy 
(density) gained by performing a phase mixture. 



The equilibrium condition (|20p for the density pro- 
file p(x) is formally equivalent to an equation of mo- 
tion for a particle of mass C moving in the potential 
V(p) — — A/(p), with p denoting the coordinate and x 
the time. (We note that Af(p) vanishes at the two co- 
existence densities and is positive in between.) Conser- 
vation of the corresponding energy ^C(d x p) 2 + V (which 



vanishes) then determines the gradient at each position 

_ 1 

2 



d x p{x) 



C 



A/(p(^)) 



(22) 



The local excess in the free energy density due to the 
interface (see App. [Cj is given by 



= Af(p(x)) + \C{d x p{x)f 



(23) 



2A/05(x)) 



The total deficit in free energy per unit interface area, 
equal to the interface tension, is then given by [23[ 



dx f To (x) 



dp(x) 



Af(p(x)) (24) 



d x p(x) 

P2 1 f 2 do 

dp[2CAf(p)} 2 = a / [2e e A f(p)Y- 

Pi Jpi Ps 



We note that this quantity can be obtained without ex- 
plicit knowledge of the interface density profile p(x) and 
it scales directly with the length parameter a. It is shown 
in Fig.[4]as a function of temperature. As expected, it de- 
creases steadily from its maximum value at T = until 
it vanishes at T c . With the (somewhat arbitrary) pa- 
rameter values adopted, the zero-temperature interface 
tension is 7g 2 « 16 MeV/fm 3 , about 16 times the famil- 
iar nuclear surface tension. This value lies near the lower 
end of the rather wide range of expected values for the 
tension between quark and nuclear matter (typical low 
values are 10— 20 MeV/fm 2 , while typical high values are 
50-100MeV/fm 2 , see for example Refs. [Il[20[). 

The density profile itself, p(x), can be obtained by in- 
tegrating Eq. 



p{x) = p{x ) + p c 



-A/(p(*)) 



- , (25) 
a 



where xq is some location where the density is known. 
We note that it would not be feasible to start the in- 
tegration at xq — > ±00, where p(xq) — > pi, since the 
gradient vanishes in the same limit, p x {xo) ^ 0. In 
stead, we take xq to be that location where the func- 
tion Af(p(x)) has its maximum. Since the derivative 
d p Af(p) = pt [p) ~ Pa thus vanishes at p = p{xq), it fol- 
lows that the local bulk chemical potential at xq matches 
the coexistence value, p?r (p(xo)) = po, and this relation 
can be used to find the starting density value p(x ). [The 
local bulk chemical potential px (p) must match the co- 
existence value pa = pT (pi) for some intermediate den- 
sity because px (p) exhibits an undulation between p\ 
and P2, going first through a maximum pt {pa) > Po 
and then through a minimum Pt (pb) < Mo, so it must 
equal po somewhere between pa and ps-} 

The density profile p(x) is shown in Fig. [5] for Tq = 
and T = \T C . At each temperature, it scales horizon- 
tally with the length parameter a. The interface profile 
can be characterized by the cumulants of the associated 
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interface location function, g(x) = p x (x)/(p2 — Pi) (see 
App. [U]). So the mean interface location is x — (x) = 
J dxxg(x), while its width b is the corresponding dis- 
persion, b 2 = ((x — x) 2 ). A convenient measure of the 
profile skewness is given by the dimensionless parame- 
ter 73 = ((x — x) 3 )/b 3 . As the temperature is increased, 
the profile grows progressively wider and more symmet- 
ric, while its mean location moves closer to x$. With 
the adopted parameter values we find x — xo = f.Ola, 
b = 2.62 a, and 73 = 0.96 at T = 0. 

It should be noted that in the present simple treat- 
ment, where the finite range is taken into account by 
means of a gradient term, the interface tension as well as 
the detailed density profile shape reflect the specific den- 
sity dependence of the free energy fr{p), i.e. they follow 
directly from the employed bulk equation of state, apart 
from scalings related to the strength of the gradient term. 

V. COLLECTIVE MODES 

We now wish to study the dynamical response to the 
introduction of small density undulations imposed on a 
system that is static and uniform, 5e(r) = i(r) — e and 
Sp(r) = p(r) — p. We first note that the local change in 
the pressure is then of a similar form, 6p(r) — p(r) — p 
with p = p(e, p) . For simplicity, we assume that the 
time evolution is described by fluid dynamics and we 
first disregard dissipation. The equations of motion then 
arise from energy-momentum conservation, 9 /J T AU/ = 0, 
together with conservation of (baryon) charge, <9 M j M = 0. 

Assuming that the local flow velocities v(r) are non- 
relativistic, we may ignore v 2 and thus put 7 to unity. 
This yields the following five equations of motion, 

= dpT^{r,t) « dtSe + hdy, (26) 
= dp,T^(r,t) » hdtrf+d^p, (27) 
= d^(r,t) a dtdp + pdy , (28) 

where h = e + p is the enthalpy density of the uniform 
system. As usual, the equations for can be combined 
to a sound-wave equation, while a comparison of the first 
and last equations yields the evolution of the density dis- 
turbance in terms of that of the energy disturbance, so 

d 2 Se(r) = d l d l 8p{r) , (29) 
hd t S P (r) = pd t Se(r) . (30) 

It is straightforward to see that, to leading order in the 
disturbances 6e(r) and Sp(r), the local pressure is 

p(r) » p(e(r),p(r))-CpV 2 p(r) . (31) 

The first term is the usual local-density approximation, 
i.e. the pressure is calculated as in uniform matter that 
has been prepared with the local density values, while 
the second term arises from the gradient correction to the 
chemical potential (fT2")> . Therefore, to the same order, 

V 2 5p(r) » p E V 2 e(r)+ Pp V 2 p(r)-CpV l p(r) , (32) 



where p e = d E p(e, p) and p p = d p p(e, p) evaluated at the 
local phase point (e, p) — (i(r),p(r)). 

If we require the undulations to be of harmonic form, 
5e(r) = Ek exp(ik-r—iujt) and Sp(r) = p^ cxp(ik-r—iujt), 
thenEq. ([30]) requires hpk — The dispersion relation 
is then readily obtained from (|2"9"|) . 

u 2 = v 2 s k 2 + ctk 4 = v 2 s k 2 + a 2 ^k A . (33) 
h h Pg 

Here the first term is what emerges in ordinary ideal fluid 
dynamics, with v s being the isentropic speed of sound 
(see Eq. (|57]l). 

v s = Ve + =V P = -^[h 2 a ee + 2hpa ep + p 2 o- pp ] , (34) 

with a sp = d e d p <r(s, p) evaluated at (e,p) = (s,p), etc.. 
This part of the dispersion relation is perfectly linear, 
LJk = v s k. That pathological behavior is modified by 
the gradient term which generally increases lu 2 . In the 
spinodal region, where v 2 is negative, the collective fre- 
quency is imaginary, ujk — ±«7fc, and the gradient term 
then suppresses the growth of high-fc modes. As a result, 
the growth rate "fk will exhibit a maximum followed by a 
rapid fall-off to zero as a function of the wave number fc, 
as is familiar from other substances exhibiting spinodal 
instability 

It is instructive to write the growth rate on the form 
7fe = |v s |fc(l — fc 2 /^max) 1 ^ 2 where the maximum wave 
number for which spinodal instability occurs is given by 

k 2 = AM = _A^IM (35) 

C p 2 £ g p 2 a 2 ' K ' 

The maximum in 7^ occurs at the "optimal" wave num- 
ber k op t = k max /y/2 and, as the amplification pro- 
cess proceeds, undulations of this size will become dom- 
inant and a characteristic spinodal pattern will thus 
emerge. The corresponding largest growth rate is 7 op t = 
^|u s |fc max = I v s I fe op t / V2- This quantity scales inversely 
with the length parameter a and the associated optimal 
wave length A op t = 27r/fc opt thus scales directly with a. 
Consequently, an increase of a will increase the scale of 
the most rapidly amplified mode as well as the associated 
shortest growth time t op t — l/7opt- 

The spinodal growth rates jk depend on the environ- 
ment, as specified for example by p and T. The tem- 
perature dependence is illustrated in Fig. [S] for p = p c , 
while the density dependence is shown in Fig.[7J 7fc(p, T) 
generally vanishes along the spinodal boundary and it de- 
creases as a function of temperature. With the present 
model, we thus find that the fastest mode at p = p c 
has a wave length of A op t ~ 3 fm and a growth time of 
^opt ~ l.Ofm/c. As the temperature is raised, the max- 
imum wave number fc max decreases as do the optimal 
values fc opt and 7 op t- While the obtained temperature 
dependence is quite significant, it should be recognized 
that the thermal properties of the present model may 
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FIG. 6: The growth rate vs. the wave number k (top) and 
the corresponding growth times tk = fo/jk vs. the wave length 
\k = 2%/k (bottom), for various temperatures T at p — po- 



FIG. 7: The optimal wave length A op t (top) for various tem- 
peratures T and the growth time tk = ft/7fc (bottom) at T = 
for several values of the wave length = 2n/k, as functions 
of the degree of compression, p/po- 



not be realistic. (By contrast, the spinodal growth rates 
in nuclear matter are relatively independent of tempera- 
ture at the low end because of the fermion nature of the 
constituents Q.) On the other hand, the dependence of 
7fc (p, T) on density is more moderate in the phase region 
of most rapid growth (as in dilute nuclear matter 0]). 

It is important to appreciate that the phase region 
of instability for ideal fluid dynamics is bounded by the 
isentropic spinodal (where v s = 0) and it therefore lies 
inside the region of thermodynamic instabilty which is 
bounded by the isothermal spinodal (see Fig. [3]). There 
are unstable isentropic modes whenever v 2 < 0. Insertion 
of the susceptibilities a ee , o~ ee , a £e (see App.[B]) into Eq. 
yields an explicit expression for the speed of sound, 



2 _ P 



d P/ S: , ,, 



IdT 



P w o 



IdT 



(36) 



so the condition for instability becomes |dT + pw'^ < 0. 
At zero temperature this amounts to dpWo(p) < 0, which 
occurs exactly within the isothermal spinodal density re- 
gion, as one would expect since T = 4=> a = 0. How- 
ever, as T is increased, the region of isentropic instabil- 
ity shrinks faster than the region of isothermal instability 
and it disappears entirely at T max = |T C . 

The above analysis was based on ideal fluid dynamics 
which conserves entropy, d^a^ = 0, where er' 1 = au^ is 
the entropy current density. We wish to conclude this 
section by briefly discussing the effects of including vis- 
cosity into the fluid-dynamic treatment. Within the non- 
relativistic framework used for the derivation of the dis- 



persion relations for the normal modes in bulk matter, 
the inclusion of shear and bulk viscosity into the fluid- 
dynamic treatment changes the pressure gradient by the 
term — V[|?7 + £]V ■ v where r\ and £ are the shear and 



bulk viscosity coefficients, respectively, 
equation is then modified accordingly, 



The dispersion 



u> 



2;„2 



vik 



h p, 



l2 £gP fc 4„ 



(37) 



where we have assumed that the combination f = |?7 + C 
can be regarded as constant, for simplicity. Clearly, the 
zero-frequency modes occur for the same wave numbers 
as before, k = and k = /c max , and the inclusion of 
viscosity does not afffect the location of the spinodal 
boundary. (However, if thermal conductivity were in- 
cluded, the spinodal boundary would gradually expand 
towards the isothermal boundary [HI].) To leading or- 
der, the viscosity adds a negative imaginary term to the 
frequency, —^k 2 /h, which in turn gives rise to an expo- 
nential damping factor. Furthermore, inside the spinodal 
region the collective frequencies are still purely imagi- 
nary, u> = iy±, and we find 



7± 



= ± 



\k 2 



• £gp_ 
e p 2 



k 4 + k 2 



h 2 



kir • (38) 



Thus the growth rate 7+ is reduced by 1=3 ^k 2 /h and 
the optimal wave number becomes smaller as well. Even 
though the qualitative features will remain the same, the 
viscous effects may be quantitatively important [13] • 
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In order to illustrate the key role played by the finite 
range in producing spinodal decomposition, let us briefly 
consider what would happen without the gradient term. 
We already noted that the resulting non-viscous disper- 
sion relation would exhibit linear growth, 7^. = l^sl^j 
and thus not favor any particular length scale. When 
viscosity is included, the growth rate would still grow 
monotonically, 9^7+ (k) > 0, but level off for large fc, 



7+(fc 



k 2 



(39) 



Thus the large-A: divergence characteristic of standard 
ideal fluid dynamics would be eliminated, but there 
would still not be a preferred length scale. 



VI. SPINODAL DECOMPOSITION? 




Density p 



Using the particular model and parameter values cho- 
sen here, we now discuss the prospects for spinodal de- 
composition to occur during a nuclear collision. Although 
the explored mean-field model is rather simplistic and the 
specific parameter values are somewhat uncertain, the re- 
sulting features appear to be within the range of plausi- 
bility. It may therefore be instructive to explore the con- 
sequences. Obviously, as further progress is made, both 
theoretically and experimentally, the models should be 
appropriately refined. 

In order to understand under what experimental con- 
ditions spinodal decomposition may actually occur, it is 
useful to consider how the thermodynamic conditions in 
the bulk of the collision system evolve in the course of 
time. Such phase trajectories were studied for gold-gold 
collisions with a variety of existing dynamical transport 
models and we shall make use of those results for 
our estimates. Ref. [l6[ calculated the evolution of the 
mechanical phase point (p(t), s(t)) in order to avoid mak- 
ing any assumption about local thermalization; since we 
are here mainly concerned with the expansion stage, we 
assume that equilibrium has been established and so we 
shall frame our discussion in terms of the canonical phase 
variables (p(t), T(t)) which are somewhat more intuitive. 

Generally speaking, the prospects for spinodal decom- 
position can be expected to be better the more time the 
bulk of the matter spends inside the region of spinodal in- 
stability. Let us therefore consider how this quantity de- 
velops with the collision energy. For the discussion below, 
we assume that the equation of state has the expected 
form with a first-order phase transition terminated by a 
critical point, as drawn schematically in Fig. [8] (see also 
Fig. [3]) . It seems natural to introduce a number of thresh- 
old values of the collision energy E: E\, Ea, Eb, E2, E c . 
Their meaning is illustrated in Fig. [8] and they will be 
explained in turn below. 

At the lowest collision energies, E < E±, the compres- 
sions achieved are insufficient to bring any part of the 
matter inside the region of phase coexistence. Conse- 



FIG. 8: Illustration of the dynamical phase traejctories for the 
most compressed matter produced at the various threshold 
collision energies Ex, Ea, Eb, E2, E c . 



quently, at such low energies, it would probably not be 
possible to probe the phase transition. 

As the collision energy is raised above E\, the phase 
trajectory (p(t),T(t)) of the most compressed matter 
makes ever larger incursions into the phase coexistence 
region. Characterizing such phase trajectories by the 
highest compression achieved, p m ax(E), we expect this 
"turning point" to gradually move across the phase co- 
existence region as E is raised. It first enters the spin- 
odal region for E = Ea and it has traversed it fully for 
E = Eb, reaching the other side of the coexistence region 
at E= E 2 . 

At collision energies above E2, the steady expansion 
of the bulk matter subsequent to its maximum compres- 
sion drags its phase trajectory through the phase coex- 
istence region (and the spinodal region within it). As E 
is increased the slope of the expansion phase trajectory 
steepens (see Ref. [16|) and the traversal time becomes 
steadily shorter, both because the expansion is faster and 
because the region of instability becomes narrower at the 
ever higher excitations encountered. At a certain "criti- 
cal" collision energy, E = E c , the phase trajectory passes 
right through the critical point (p c ,T c ) and at supercrit- 
ical collision energies, E > E c , the phase trajectory will 
miss the unstable phase region altogether. 

Generally, the evolving local thermodynamic condi- 
tions during a collision will differ from one location to 
another. Consequently, a single collision event gives rise 
to an entire bundle of phase trajectories and the above 
discussion pertains to just the phase trajectory of the 
most compressed matter of the collision system which, 
for a symmetric collision, is presumably located around 
the center. Furthermore, there is a dependence on the 
geometric features of the collision system, such as the nu- 
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clear sizes and the impact parameter. Thus the precise 
meaning of the various threshold energies is somewhat 
fuzzy and they play primarily a conceptual role. This 
underscores the fact that quantitative predictions must 
rely on detailed dynamical calculations. 

Our special interest here concerns the relatively nar- 
row interval of collision energy within which the turning 
point lies inside the spinodal region, Ea < E < Eb- In- 
tuitively, one would expect that collision energies slightly 
below Eb would be optimal for maximizing the time 
spent by the phase trajectory inside the spindoal region. 
Such collisions, in turn, would presumably be most favor- 
able for the development of spinodal decompostion. The 
transport calculations reported in Ref. [l6| suggest that 
this optimal beam energy is 5 — 15 GeV per nucleon for 
a stationary target setup. The presice values depend not 
only on the specific location of the spinodal phase bound- 
aries, i.e. on the specific equation of state (which is still 
unknown) , but also on the complications arising from the 
non-uiformity of the density and its time evolution. 

In order to get a rough idea of the degree of spin- 
odal growth that could be expected during a given colli- 
sion, we assume that we know the time evolution of the 
bulk density, p(t), and the associated temperature, T(t). 
As pointed out above, these quantities are local and the 
present analysis employs suitable average values that can 
be taken as representative of an extended part of the sys- 
tem. For a given collision energy between Ea and Eb, 
the phase trajectory (p(t),T(t)) enters the spinodal re- 
gion at the time t = t{ and exits it again at the time 
t = tf. Thus, for ti < t < tf the collective dispersion re- 
lation yields unstable modes with associated growth rates 
7fc(i) = 7fc(p(t),T(t)) which serve to amplify irregulari- 
ties in the density. 

An accurate calculation would need to take account 
not only of the distribution of density fluctuations but 
also of the fact the entire scenario changes in time (and 
relatively rapidly) . This is beyond our present scope and 
we seek to obtain a simple estimate by considering the 
following amplification coefficient [Til ]. 

To = J% Q {t)dt « §7maxAt, (40) 

where 70 (t) is the maximum growth rate at the time t, 
7o(*) = lo P t{p(t),T(t)) and 7 max is the largest growth 
rate overall. The factor of two thirds accounts roughly 
for the fact that 70 (t) has a parabola- like appearance, 
starting out from zero at U, exhibiting a broad maximum 
of 7max, and then dropping to zero again at tt, so we put 
(70 (i)) ~ §7max- Judging from the transport calculations 
reported in Ref. [IB], we estimate the duration of the 
spinodal stage to be At = tf — U w 6fm/c. To estimate 
7max is more difficult. Our present calculations give an 
overall fastest growth time of » lfm/c, obtained for rel- 
atively broad range of densties and for zero temperature. 
However, the compression achieved in a nuclear collision 
is inevitably accompanied by a corresponding agitation, 



and we therefore expect T/T c = | — i to be more re- 
alistic. In this connection it should be realized that the 
the rapid decrease of the calculated growth rate as the 
temperature is increased is to some extent a reflection of 
the fact that the present instability region is bounded the 
isentropic rather than the isothermal spinodal line. On 
the other hand, the inclusion of dissipation (which would 
expand the region of instability) is expected to slow the 
dynamics down to a significant degree [l7| . Therefore it 
is probably more realistic to expect the fastest growth 
time to be several times that most optimistic value, so 
we use 7^ x « 2 -4fm/c 

With these rough numbers, we then find the value of 
the amplification coefficient to be Tq sa |-6/ (2—4) = 1 — 2. 
The corresponding amplitude growth factor [14| is then 
given by Go = exp(r ) ~ 2.7— 7 .4. When trying to judge 
the significance of this value, one should keep in mind 
that the density-density correlation function is propor- 
tional to the square of the amplitude growth factor, i.e. 
(5p(n) exp(iri2-k)8p(r 2 )) ~ G\. 

It should also be realized that after the system ex- 
its the spinodal instability region of the phase diagram, 
it has still to traverse the metastable region between 
the spinodal boundary and the phase coexistence line. 
While nearly uniform matter is mechanically stable in 
this regime, this is no longer so for matter having signif- 
icant deviations from uniformity. Therefore, the undula- 
tions resulting from even relatively modest amplifications 
during the spinodal stage may be further amplified dur- 
ing the metastable stage and thus lead to observationally 
interesting clumping of the expanding matter. 

The above numerical estimates were obtained for the 
adopted range value of a = 0.2 fm which is of course 
rather uncertain. We recall that it leads to an interface 
tension of 7y= ~ 16MeV/fm 2 and an optimal wave- 
length of A op t ~ 4fm at T ~ \T C . The interface tension 
is at the lower end of what has been used by various 
authors [H, H3], but does not appear to be unreason- 
able considering the large uncertainties on this quantity. 
As for the wave length, it is of interest to note that a 
spherical volume having such a diameter would contain 
a baryon number of Bq = fAgPc w 24. At a temperature 
of half the critical value, the completion of the phase de- 
composition would distribute this matter approximately 
evenly between the two coexisting phases, assuming the 
high-density phase is concentrated in a sphere embedded 
into the low-density phase. Such a blob of deconfmcd 
matter is large enough to constitute a macroscopic sta- 
tistical source, while at the same time being probably a 
sufficiently small part of the total system to permit the 
simultaneous formation of several such blobs and thus 
make it feasible to perform a size-correlation analysis. 

If we were to use only half that range, a = 0.1 fm, those 
quantities would decrease correspondingly to jj?=o ~ 
8MeV/fm , which would be somewhat low in compar- 
ison to the existing estimates, though perhaps not im- 
possible, and A opt (^T c ) w 2fm, leading to Bo w 3 which 
seems too small to constitute a macroscopic source that 
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could have observational significance. On the other hand, 
if we were to double the range, a — 0.4fm/c, we would 
obtain 7j?=o ~ 32 MeV/fm 2 , which would seem quite rea- 
sonable, but A op t(|T c ) w 8fm (hence Bq ss 192) would 
be far too large to produce a useful effect. 

This analysis reveals that whether a given model leads 
to spinodal phase decomposition in simulations of nuclear 
collisions depends rather delicately on the specific values 
of its parameters. Consequently, it is far from certain 
that collisions of real nuclei would in fact produce this 
phenomenon at any collision energy. 

The estimates above suggest that collisions within a 
suitably tuned energy range may produce bulk matter 
that stays inside the mechanically unstable phase re- 
gion sufficiently long for some degree of spinodal clump- 
ing to occur. However, it is hard to predict whether 
the amplification of the fluctuations will suffice to bring 
about the characteristic spinodal enhancement of a cer- 
tain length scale. Indeed, the degree of amplification 
might be so marginal that the matter will revert to ap- 
proximate uniformity after reentering the stable regime 
and no clumping- like phase separation would then occur. 

This uncertainty underscores the need for studying the 
observable consequences of a spinodal decomposition. In- 
deed, only if the spinodal phenomenon manifests itself in 
detectable signals can it be turned into a useful tool for 
probing the equation of state. It is as of yet far from clear 
whether any proposed "signals" of the phase decomposi- 
tion would in fact survive the subsequent expansion stage 
dominated by hadronic resonances. While current mis- 
croscopic transport models present useful tools for such 
investigations, the insight that can be gained will likely 
remain somewhat limited until we achieve a better un- 
derstanding of the phase transition dynamics itself. 



VII. CONCLUDING REMARKS 

The present study was motivated by the need for devel- 
oping theoretical models that can address the dynamics 
of the confinement phase transition that is expected to 
occur during the expansion stage of collisions between 
heavy nuclei at suitably tuned energies. As a step to- 
wards this goal, we have developed a simple model within 
which we have studied both the collective dispersion rela- 
tion for the mechanically unstable modes of bulk matter 
in the spinodal region of the thermodynamic phase dia- 
gram and the properties of the interface between the two 
coexisting phases into which such an ustable system seeks 
to decompose. These properties are central to the phase 
transition dynamics and since they are fundamentally re- 
lated it is important they that be treated consistently. As 
far as we are aware, this is the first time that these dif- 
ferent properties have been addressed within the same 
model framework. 

The key element is the inclusion of a finite interaction 
range without which there would be neither an inter- 
face tension nor spinodal decomposition. Indeed, a zero- 



range model would render any interface perfectly sharp 
and there would be no associated energy cost, hence no 
basis for determining the geometric structure of a system 
composed of coexisting phases (a zero-range model would 
admit even fractal intermingling of the phases). Further- 
more, without a finite range to suppress the dynamics 
of short wavelength disturbances, the growth rate would 
increase steadily with wave number (even in the presence 
of viscosity) , hence not dispay a maximum as is charac- 
teristic of spinodal decomposition. Thus, any model that 
aspires to be of use for phase transition dynamics must 
incorporate a finite range. 

A recent study by Skokov and Voskresensky has sought 
to accomplish that by means of a gradient term in the 
free energy functional [TtJ • In the present study, we take 
a similar approach by constructing an equation of state 
in which the interaction-energy density in non-uniform 
matter contains a gradient term. This can be thought 
of as approximating a convolution with a kernel of finite 
range, as is often done for the mean field in low-energy 
nuclear physics. 

Considering that no quantitatively reliable calculations 
are yet available, we have had to fix the parameters of the 
equation of state for bulk matter on the basis of our best 
guess for the location of the coexistence region, includ- 
ing the critical point. Consequently, our results should 
not be considered as more than suggestive. The addi- 
tional range parameter entering into the gradient term 
has been adjusted to yield reasonable values for the inter- 
face tension and the spinodal growth rates. Importantly, 
a change in this range by a factor of two or more would 
render the model results either implausible (relative to 
existing estimates) or phenomenologically uninteresting 
(in that the resulting model would not produce spinodal 
decomposition in a collision scenario). 

With the gradient term included, we have then studied 
the equilibrium interface between two coexisting phases 
and determined the temperature dependence of the den- 
sity profile and the associated interface tension. 

In order to address the collective modes in bulk matter, 
a dynamcial model is needed and we have adopted fluid 
dynamics, which can readily be adapted to the finite- 
range equation of state. The gradient term suppresses the 
growth of short wavelengths and thus yields a physically 
reasonable dispersion relation for the spinodal modes. 

Taking guidance from existing phase trajectories ex- 
tracted from various transport simulations [16[, we have 
used the calculated growth rates to estimate the degree 
of amplification that might occur when the collision en- 
ergy is adjusted to maximize the exposure to the spinodal 
instabilities. The resulting amplification amounts to one 
or two factors of e, which may suffice to trigger a phase 
separation due to the subsequent further amplification 
from the intermediate metastable phase region. While 
this conclusion gives grounds for guarded optimism, it 
also brings out the fact that a full dynamical simulation 
is needed for a more detailed assessment. 
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APPENDIX A: COMPRESSIONAL ENERGY 



# 


T (MeV) p/po 


e(MeV/fm 3 ) 


p (MeV/ fin 3 ) 


p (MeV) 


1 





3.00 


182 


171 


769 


A 





3.76 


276 


201 


829 


c 


170 


4.70 


729 


395 


923 


B 





6.18 


563 


130 


732 


2 





8.01 


772 


171 


769 



TABLE I: Values of the temperature T, compression p/po, 
energy density e, pressure p, and chemical potential p at the 
two zero-temperature coexistence points (#1 and #2), the 
two spinodal boundaries at zero temperature (#A and #B), 
and the critical point (#c) (see Figs. [T] and [2}. 



The key quantity in the employed illustrative equation 
of state is the energy density associated with the com- 
pression of bulk matter at zero temperature, w (p) = 
£t=o{p)- We obtain this basic function by interpolating 
between a "hadron gas" and a "quark-gluon plasma" , 

w o{p) = x(p)w H {p) + [1 -x(p)]wq(p) , (Al) 

where the interpolation function is taken as 

X (p) = [l + e^-P^/^}- 1 (A2) 

and each phase is taken to display a simple power form, 

w H (p) = c H (p/p ) 2 , w Q (p) = c Q (p/p ) 4/3 + S, (A3) 

Po ~ 0.153 fm -3 being the nuclear saturation density. 

For a wide range of parameter values, the resulting 
compressional energy density wo(p) exhibits a region of 
negative curvature, thus ensuring the existence of a first- 
order phase transition. [We recall that phase coexistence 
requires equal chemical potentials, pi(pi) = p{ P 2), hence 
equal slopes of u>o(p) (since pr=o = dpWa), as well as 
equal pressures, p(pi) = p(p2), hence a common tangent 
of w (p) (since p T =o = PP~ f = pd p w - w ).] Adopting 
the specific values ch — 92.6, cq = 288.9, and B = 
408.3 (all in McV/fm 3 ) together with p x = 3.2535 p and 
p w = 0.945 poj w e obtain p\ = 3/?o and P2 = &Po for 
the zero-temperature coexistence densities. (It would be 
straightforward to obtain other values by readjusting the 
parameters.) 

The resulting compressional energy density Wo(p) is 
shown in Fig. Q] by the solid (black) curve, while the 
dashed (blue and red) curves are the individual functions 
wh(p) and wq(p). It is apparent from the plot that the 
coexistence features depend delicately on the parameter 
values. One would therefore expect any model calcula- 
tion of the phase structure to be endowed with rather 
large uncertainties, thus reinforcing the need for experi- 
mental information. 

Table |T] summarizes the values of the various quantities 
of interest at the phase and spinodal boundaries: 



APPENDIX B: SUSCEPTIBILITIES 

Mechanical stability is determined by the curvature 
tensor cr of the entropy density, which in the present 
model has the elements cr e6 , cr pe — a £p , a pp , where 



dMe, P ) 



P 



d p a(s,p) 



2 (e - w ) 2 - 1; 

1 dP^O 
l< 2 K 

. d ,w{ 

K 



14 < °» ( B1 ) 



a p£ = d p d s a(£, P) = °e P = lz + 1^ > (B2) 



P 



d ro_4^o_ |^L.(B3) 



K 



with k(e, p) = e — wn(p) = \dpT. We generally have 



pT 



dpPT(p) 



which in the present case amounts to 
2 \dT + pw'J 



(B4) 



(B5) 



Since the occurrence of mechanical instability requires 
that at least one of the eigenvalues of a be positive, and 
| cr | is the product of the eigenvalues, it follows that me- 
chanical instability occurs at densities for which w'q(p) 
is negative and then extends up to the temperature 
T max (p) = — 2P w aip)- The critical temperature is the 
largest of those, T c = T max (p c ) = ^p c w'^(p c ). 
The isothermal sound speed vt readily follows, 



p_fdp s 
h \dpJ T 



T,» T - 

h <j ff 



JdT^r . (B6) 



where h = p + e and we have used (dp/dp)r — dpPr{p), 
while the isentropic sound speed v s is given by 



[h 2 a ee + 2hpa ep + p 2 a pp ] 



(B7) 
dT + pw'^ 



IdT 



where s = a/p is the entropy per particle and we have 
used that the requirement Ss = implies pSe = h5p. We 
note that v 2 s > v\ for T > 0. 
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APPENDIX C: INTERFACE 

Consider a planar interface between two semi-infinite 
systems and let the bulk values of the charge and energy 
densities in system i be and pi, respectively, and as- 
sume that p\ < pi- Presumably the local densities p(x) 
and e(x) rapidly approach these asymptotic values away 
from the interface. We may generally define the interface 
location function [2l[ , 

g p (x) = , [ °°g P (x)dx = 1 , (CI) 

which peaks near the interface and is normalized to unity. 
Its moments provide quantitative characteristics of the 
interface profile and we list here the first three [2l| . 

/+00 
xg p (x)dx , (C2) 
-00 

Width: b = [((x - xf)] 1/2 , (C3) 
Skewness: 73 = ((x-x) 3 )/b 3 . (C4) 



In order to extract the interface tension, we follow the 
discussion in Ref. [23]. Thus, the difference between the 
actual diffuse energy density e(x) and that associated 
with the corresponding sharp-surface configuration that 
would result if there were no gradient term is given by 

£12(2;) = e{x)-Ei {p(x)-pi), 1 = 1,2. (C5) 

92 - pi 

The function £12(2;) is obviously smooth and, moreover, 
it tends quickly to zero away from the interface, so it is 
indeed peaked in the interface region. (It is elementary 
to see that it does not matter whether one uses i = 1 or 
i = 2 in the above expression.) Its integral is then the 
total energy deficit associated with the diffuse interface, 
the interface tension, 

/+00 
£12(2:) dx , (C6) 
-00 

which can thus readily be calculated once the profile den- 
sities p(x) and e(x) are known. 
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